#include <iostream>
#include <opencv2/core.hpp>
#include <opencv2/imgcodecs.hpp>
#include <opencv2/imgproc.hpp>
#include <opencv2/highgui.hpp>

static void _get_a_b_c(int distType, int maskSize, float& a, float& b, float& c) {
    switch (distType) {
        case cv::DIST_C:
            CV_Assert(maskSize == cv::DIST_MASK_3);
            a = 1.0f;
            b = 1.0f;
            break;
        case cv::DIST_L1:
            CV_Assert(maskSize == cv::DIST_MASK_3);
            a = 1.0f;
            b = 2.0f;
            break;
        case cv::DIST_L2:
            if (maskSize == cv::DIST_MASK_3) {
                a = 0.955f;
                b = 1.3693f;
            } else if (maskSize == cv::DIST_MASK_5) {
                a = 1.f;
                b = 1.4f;
                c = 2.1969f;
            }
            break;
        default:
            CV_Error(cv::Error::StsBadArg, "Invalid Distance Type!");
    }
}

static const int DIST_SHIFT = 16;
static const int INIT_DIST = INT_MAX;
static const int MAX_DIST = (INT_MAX >> 2);

static unsigned int _float_2_fix_unsigned_int(float x, int shift) {
    return cvRound(x * (1 << shift));
}

static void _init_border(cv::Mat& mat, int borderSize) {
    int rows = mat.rows, cols = mat.cols;
    for (int i = 0 ; i < borderSize; i++) {
        int* pTop = mat.ptr<int>(i);
        int* pBottom = mat.ptr<int>(rows-1-i);
        for (int j = 0; j < cols; j++) {
            pTop[j] = INIT_DIST;
            pBottom[j] = INIT_DIST;
            if (j < borderSize || j > cols - borderSize-1) {
                mat.ptr<int>(i)[j] = INIT_DIST;
            }
        }
    }
}

static cv::Mat _distance_transform_3x3(const cv::Mat& src, int distType) {
    const int BORDER_SIZE = 1;
    int maskSize = cv::DIST_MASK_3;
    float a, b, c;
    _get_a_b_c(distType, maskSize, a, b, c);
    
    const unsigned int HV_DIST = _float_2_fix_unsigned_int(a, DIST_SHIFT);
    const unsigned int DIAG_DIST = _float_2_fix_unsigned_int(b, DIST_SHIFT);
    const float SCALE =1.f / (1 << DIST_SHIFT);
    
    cv::Mat temp(src.rows+BORDER_SIZE*2, src.cols+BORDER_SIZE*2, CV_32SC1);
    _init_border(temp, BORDER_SIZE);
    
    cv::Mat dist(src.size(), CV_32F);
    
    // forward pass
    for (int i = 0; i < src.rows; i++) {
        for (int j = 0; j < src.cols; j++) {
            int tempRowCenter = i + BORDER_SIZE;
            int tempColCenter = j + BORDER_SIZE;
            if (src.ptr<uchar>(i)[j] == 0) {
                temp.ptr<int>(tempRowCenter)[tempColCenter] = 0;
            } else {
                unsigned int t0 = temp.ptr<int>(tempRowCenter-1)[tempColCenter-1] + DIAG_DIST;
                unsigned int t = temp.ptr<int>(tempRowCenter-1)[tempColCenter] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter)[tempColCenter-1] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                temp.ptr<int>(tempRowCenter)[tempColCenter] = t0;
            }
        }
    }
    
    // backward pass
    for (int i = src.rows-1; i >= 0; i--) {
        for (int j = src.cols-1; j >= 0; j--) {
            int tempRowCenter = i + BORDER_SIZE;
            int tempColCenter = j + BORDER_SIZE;
            unsigned int t0 = temp.ptr<int>(tempRowCenter)[tempColCenter];
            if (t0 > HV_DIST) {
                unsigned int t = temp.ptr<int>(tempRowCenter+1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter-1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter)[tempColCenter+1] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                temp.ptr<int>(tempRowCenter)[tempColCenter] = t0;
            }
            t0 = (t0 > MAX_DIST ? MAX_DIST : t0);
            dist.ptr<float>(i)[j] = (float)(t0 * SCALE);
        }
    }
    return dist;
}

static cv::Mat _distance_transform_for_L1_8U(const cv::Mat& src) {
    CV_Assert(src.type() == CV_8UC1);
    
    cv::Mat dist(src.size(), CV_8U);
    
    uchar lut[256];
    for (int i = 0; i < 256; i++) {
        lut[i] = cv::saturate_cast<uchar>(i+1);
    }
    
    // forward scan
    dist.ptr<uchar>(0)[0] = (uchar)(src.ptr<uchar>(0)[0] == 0 ? 0 : 255);
    
    for (int i = 1; i < src.cols; i++) {
        dist.ptr<uchar>(0)[i] = (uchar)(src.ptr<uchar>(0)[i] == 0 ? 0 : lut[dist.ptr<uchar>(0)[i-1]]);
    }
    
    for (int i = 1; i < src.rows; i++) {
        int a = (src.ptr<uchar>(i)[0] == 0 ? 0 : lut[dist.ptr<uchar>(i-1)[0]]);
        dist.ptr<uchar>(i)[0] = (uchar)a;
        for (int j = 1; j < src.cols; j++) {
            a = (src.ptr<uchar>(i)[j] == 0 ? 0 : lut[MIN(a, dist.ptr<uchar>(i-1)[j])]);
            dist.ptr<uchar>(i)[j] = (uchar)a;
        }
    }
    
    // backward scan
    int a = dist.ptr<uchar>(src.rows-1)[src.cols-1];
    for (int i = src.cols-2; i >= 0; i--) {
        a = lut[a];
        dist.ptr<uchar>(src.rows-1)[i] = (uchar)(MIN((uchar)a, dist.ptr<uchar>(src.rows-1)[i]));
    }
    for (int i = src.rows-2; i >= 0; i--) {
        a = lut[dist.ptr<uchar>(i+1)[src.cols-1]];
        a = dist.ptr<uchar>(i)[src.cols-1] = (uchar)(MIN(a, dist.ptr<uchar>(i)[src.cols-1]));
        for (int j = src.cols-2; j >= 0; j--) {
            int b = dist.ptr<uchar>(i+1)[j];
            a = lut[MIN(a, b)];
            a = MIN(a, dist.ptr<uchar>(i)[j]);
            dist.ptr<uchar>(i)[j] = (uchar)a;
        }
    }
    return dist;
}

static cv::Mat _distance_transform_5x5(const cv::Mat& src, int distType) {
    const int BORDER_SIZE = 2;
    int maskSize = cv::DIST_MASK_5;
    float a, b, c;
    _get_a_b_c(distType, maskSize, a, b, c);
    
    const unsigned int HV_DIST = _float_2_fix_unsigned_int(a, DIST_SHIFT);
    const unsigned int DIAG_DIST = _float_2_fix_unsigned_int(b, DIST_SHIFT);
    const unsigned int LONG_DIST = _float_2_fix_unsigned_int(c, DIST_SHIFT);
    const float SCALE =1.f / (1 << DIST_SHIFT);
    
    cv::Mat temp(src.rows+BORDER_SIZE*2, src.cols+BORDER_SIZE*2, CV_32SC1);
    _init_border(temp, BORDER_SIZE);
    
    cv::Mat dist(src.size(), CV_32F);
    
    // forward pass
    for (int i = 0; i < src.rows; i++) {
        for (int j = 0; j < src.cols; j++) {
            int tempRowCenter = i + BORDER_SIZE;
            int tempColCenter = j + BORDER_SIZE;
            if (src.ptr<uchar>(i)[j] == 0) {
                temp.ptr<int>(tempRowCenter)[tempColCenter] = 0;
            } else {
                unsigned int t0 = temp.ptr<int>(tempRowCenter-2)[tempColCenter-1] + LONG_DIST;
                unsigned int t = temp.ptr<int>(tempRowCenter-2)[tempColCenter+1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter)[tempColCenter-1] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                temp.ptr<int>(tempRowCenter)[tempColCenter] = t0;
            }
        }
    }
    
    // backward pass
    for (int i = src.rows-1; i >= 0; i--) {
        for (int j = src.cols-1; j >= 0; j--) {
            int tempRowCenter = i + BORDER_SIZE;
            int tempColCenter = j + BORDER_SIZE;
            unsigned int t0 = temp.ptr<int>(tempRowCenter)[tempColCenter];
            if (t0 > HV_DIST) {
                unsigned int t = temp.ptr<int>(tempRowCenter+2)[tempColCenter+1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+2)[tempColCenter-1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter+2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter-1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter-2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                t = temp.ptr<int>(tempRowCenter)[tempColCenter+1] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                }
                temp.ptr<int>(tempRowCenter)[tempColCenter] = t0;
            }
            t0 = (t0 > MAX_DIST ? MAX_DIST : t0);
            dist.ptr<float>(i)[j] = (float)(t0 * SCALE);
        }
    }
    return dist;
}

static cv::Mat _true_distance_transform(const cv::Mat& src) {
    CV_Assert(src.type() == CV_8UC1);
    const float inf = 1e15f;
    int rows = src.rows, cols = src.cols;
    
    cv::Mat dist(src.size(), CV_32F);
    
    std::vector<int> drow(rows);
    for (int col = 0; col < cols; col++) {
        int distance = rows-1;
        for (int j = 0; j < rows; j++) {
            if (src.ptr<uchar>(j)[col] == 0) {
                distance = 0;
            } else {
                distance += 1;
            }
            drow[j] = distance;
        }
        distance = rows-1;
        for (int j = rows-1; j >= 0; j--) {
            distance += 1;
            if (distance > drow[j]) {
                distance = drow[j];
            }
            dist.ptr<float>(j)[col] = distance * distance;
        }
    }
    
    std::vector<int> v(cols);
    std::vector<float> z(cols+1), f(cols);
    for (int row = 0; row < rows; row++) {
        int k = 0, p;
        v[0] = 0;
        z[0] = -inf;
        z[1] = inf;
        f[0] = dist.ptr<float>(row)[0];
        for (int q = 1; q < cols; q++) {
            f[q] = dist.ptr<float>(row)[q];
            float fq = f[q];
            while (true) {
                p = v[k];
                float fvk = dist.ptr<float>(row)[p];
                float s = (fq + q*q - fvk - p*p) / (2*(q-p));
                if (s > z[k]) {
                    k++;
                    v[k] = q;
                    z[k] = s;
                    z[k+1] = inf;
                    break;
                }
                k--;
            }
        }
        k = 0;
        for (int q = 0; q < cols; q++) {
            while (z[k+1] < q) {
                k++;
            }
            p = v[k];
            dist.ptr<float>(row)[q] = std::sqrt((q-p)*(q-p) + f[p]);
        }
    }
    return dist;
}

static cv::Mat _distance_transform_5x5_ex(const cv::Mat& src, int labelType, cv::Mat& labels) {
    CV_Assert(src.type() == CV_8UC1);
    CV_Assert(labelType == cv::DIST_LABEL_CCOMP || labelType == cv::DIST_LABEL_PIXEL);
    
    labels = cv::Mat(src.size(), CV_32S);
    labels.setTo(cv::Scalar::all(0));
    if (labelType == cv::DIST_LABEL_CCOMP) {
        cv::Mat zeroPix = (src == 0);
        cv::connectedComponents(zeroPix, labels, 8, CV_32S, cv::CCL_WU);
    } else {
        int k = 1;
        for (int i = 0; i < src.rows; i++) {
            for (int j = 0; j < src.cols; j++) {
                if (src.ptr<uchar>(i)[j] == 0) {
                    labels.ptr<int>(i)[j] = k;
                    k++;
                }
            }
        }
    }
    
    const int BORDER_SIZE = 2;
    int maskSize = cv::DIST_MASK_5;
    float a, b, c;
    _get_a_b_c(cv::DIST_L2, maskSize, a, b, c);
    
    const unsigned int HV_DIST = _float_2_fix_unsigned_int(a, DIST_SHIFT);
    const unsigned int DIAG_DIST = _float_2_fix_unsigned_int(b, DIST_SHIFT);
    const unsigned int LONG_DIST = _float_2_fix_unsigned_int(c, DIST_SHIFT);
    const float SCALE =1.f / (1 << DIST_SHIFT);
    
    cv::Mat temp(src.rows+BORDER_SIZE*2, src.cols+BORDER_SIZE*2, CV_32SC1);
    _init_border(temp, BORDER_SIZE);
    
    cv::Mat dist(src.size(), CV_32F);
    
    // forward pass
    for (int i = 0; i < src.rows; i++) {
        for (int j = 0; j < src.cols; j++) {
            int tempRowCenter = i + BORDER_SIZE;
            int tempColCenter = j + BORDER_SIZE;
            if (src.ptr<uchar>(i)[j] == 0) {
                temp.ptr<int>(tempRowCenter)[tempColCenter] = 0;
            } else {
                unsigned int t0 = INIT_DIST, t;
                int label = 0;
                t = temp.ptr<int>(tempRowCenter-2)[tempColCenter-1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-2)[j-1];
                }
                t = temp.ptr<int>(tempRowCenter-2)[tempColCenter+1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-2)[j+1];
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter-2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-1)[j-2];
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter-1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-1)[j-1];
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-1)[j];
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-1)[j+1];
                }
                t = temp.ptr<int>(tempRowCenter-1)[tempColCenter+2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i-1)[j+2];
                }
                t = temp.ptr<int>(tempRowCenter)[tempColCenter-1] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i)[j-1];
                }
                temp.ptr<int>(tempRowCenter)[tempColCenter] = t0;
                labels.ptr<int>(i)[j] = label;
            }
        }
    }
    
    // backward pass
    for (int i = src.rows-1; i >= 0; i--) {
        for (int j = src.cols-1; j >= 0; j--) {
            int tempRowCenter = i + BORDER_SIZE;
            int tempColCenter = j + BORDER_SIZE;
            unsigned int t0 = temp.ptr<int>(tempRowCenter)[tempColCenter];
            int label = labels.ptr<int>(i)[j];
            if (t0 > HV_DIST) {
                unsigned int t = temp.ptr<int>(tempRowCenter+2)[tempColCenter+1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+2)[j+1];
                }
                t = temp.ptr<int>(tempRowCenter+2)[tempColCenter-1] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+2)[j-1];
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter+2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+1)[j+2];
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter+1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+1)[j+1];
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+1)[j];
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter-1] + DIAG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+1)[j-1];
                }
                t = temp.ptr<int>(tempRowCenter+1)[tempColCenter-2] + LONG_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i+1)[j-2];
                }
                t = temp.ptr<int>(tempRowCenter)[tempColCenter+1] + HV_DIST;
                if (t0 > t) {
                    t0 = t;
                    label = labels.ptr<int>(i)[j+1];
                }
                temp.ptr<int>(tempRowCenter)[tempColCenter] = t0;
                labels.ptr<int>(i)[j] = label;
            }
            t0 = (t0 > MAX_DIST ? MAX_DIST : t0);
            dist.ptr<float>(i)[j] = (float)(t0 * SCALE);
        }
    }
    return dist;
}

int main(int argc, char **argv) {
    cv::Mat src = cv::imread("/home/xlll/Downloads/opencv/samples/data/cards.png", cv::IMREAD_COLOR);
    if (src.empty()) {
        std::cout << "failed to read image" << std::endl;
        return EXIT_FAILURE;
    }
    for (int i = 0; i < src.rows; i++) {
        for (int j = 0; j < src.cols; j++) {
            if (src.ptr<cv::Vec3b>(i)[j] == cv::Vec3b(255, 255, 255)) {
                src.ptr<cv::Vec3b>(i)[j] = cv::Vec3b(0, 0, 0);
            }
        }
    }
    
    cv::Mat gray;
    cv::cvtColor(src, gray, cv::COLOR_BGR2GRAY);
    cv::Mat bw;
    cv::threshold(gray, bw, 0, 255, cv::THRESH_BINARY | cv::THRESH_OTSU);
    
    /*
    cv::Mat dist1, dist2;
    cv::Mat labels1, labels2;
    //cv::distanceTransform(bw, dist1, labels1, cv::DIST_L2, cv::DIST_MASK_5, cv::DIST_LABEL_CCOMP);
    cv::distanceTransform(bw, dist1, labels1, cv::DIST_L2, cv::DIST_MASK_5, cv::DIST_LABEL_PIXEL);
    //dist2 = _distance_transform_5x5_ex(bw, cv::DIST_LABEL_CCOMP, labels2);
    dist2 = _distance_transform_5x5_ex(bw, cv::DIST_LABEL_PIXEL, labels2);
    cv::Mat diff = dist1 - dist2;
    cv::Mat diff2 = labels1 - labels2;
    std::vector<cv::Point> nonzeros, nonzeros2;
    cv::findNonZero(diff, nonzeros);
    cv::findNonZero(diff2, nonzeros2);
    std::cout << "Different pixel number for L2 5x5: " << nonzeros.size() << std::endl;
    std::cout << "Different label number for L2 5x5: " << nonzeros2.size() << std::endl;
    */
    /*
    cv::Mat dist1, dist2;
    cv::distanceTransform(bw, dist1, cv::DIST_L2, cv::DIST_MASK_PRECISE, CV_32F);
    dist2 = _true_distance_transform(bw);
    cv::Mat diff = dist1 - dist2;
    std::vector<cv::Point> nonzeros;
    cv::findNonZero(diff, nonzeros);
    std::cout << "Different pixel number for true distance transform: " << nonzeros.size() << std::endl;
    for (int i = 0; i < nonzeros.size(); i++) {
        std::cout << nonzeros[i] << "\t";
        std::cout << dist1.ptr<float>(nonzeros[i].y)[nonzeros[i].x] << "\t";
        std::cout << dist2.ptr<float>(nonzeros[i].y)[nonzeros[i].x] << "\t";
        std::cout << diff.ptr<float>(nonzeros[i].y)[nonzeros[i].x] << std::endl;
    }
    */
    
    cv::Mat dist1, dist2;
    cv::distanceTransform(bw, dist1, cv::DIST_L2, cv::DIST_MASK_5, CV_32F);
    dist2 = _distance_transform_5x5(bw, cv::DIST_L2);
    cv::Mat diff = dist1 - dist2;
    std::vector<cv::Point> nonzeros;
    cv::findNonZero(diff, nonzeros);
    std::cout << "Different pixel number for L2 (5x5): " << nonzeros.size() << std::endl;
    for (int i = 0; i < nonzeros.size(); i++) {
        std::cout << nonzeros[i] << "\t";
        std::cout << dist1.ptr<float>(nonzeros[i].y)[nonzeros[i].x] << "\t";
        std::cout << dist2.ptr<float>(nonzeros[i].y)[nonzeros[i].x] << "\t";
        std::cout << diff.ptr<float>(nonzeros[i].y)[nonzeros[i].x] << std::endl;
    }
    
    /*
    cv::Mat dist1, dist2;
    cv::distanceTransform(bw, dist1, cv::DIST_L1, cv::DIST_MASK_3, CV_8U);
    dist2 = _distance_transform_for_L1_8U(bw);
    cv::Mat diff = dist1 - dist2;
    std::vector<cv::Point> nonzeros;
    cv::findNonZero(diff, nonzeros);
    std::cout << "Different pixel number for L1 8U: " << nonzeros.size() << std::endl;
    */
    /*
    std::vector<int> distTypes = {cv::DIST_C, cv::DIST_L1, cv::DIST_L2};
    std::vector<std::string> distTypeNames = {"C(Chessboard)", "L1(City Block)", "L2(Euclidean)"};
    std::vector<cv::Mat> dists1, dists2;
    for (unsigned int i = 0; i < distTypes.size(); i++) {
        cv::Mat dist1, dist2;
        cv::distanceTransform(bw, dist1, distTypes[i], cv::DIST_MASK_3, CV_32F);
        dists1.push_back(dist1);
        dist2 = _distance_transform_3x3(bw, distTypes[i]);
        dists2.push_back(dist2);
    }
    
    std::vector<cv::Mat> diffs;
    for (unsigned i = 0; i < dists1.size(); i++) {
        diffs.push_back(dists1[i] - dists2[i]);
    }
    for (unsigned int i = 0; i < diffs.size(); i++) {
        std::vector<cv::Point> nonzeros;
        cv::findNonZero(diffs[i], nonzeros);
        std::cout << "Different pixel number for " << distTypeNames[i] << ": " << nonzeros.size() << std::endl;
    }
    */
    /*
    std::vector<int> maskSizes = {cv::DIST_MASK_3, cv::DIST_MASK_5, cv::DIST_MASK_PRECISE};
    std::vector<std::string> winNames = {"3", "5", "PRECISE"};
    std::vector<cv::Mat> dsts, normDsts;
    for (int i = 0; i < maskSizes.size(); i++) {
        cv::Mat dst, normDst;
        //cv::distanceTransform(bw, dst, cv::DIST_L1, maskSizes[i], CV_32F);
        //cv::distanceTransform(bw, dst, cv::DIST_L2, maskSizes[i], CV_32F);
        cv::distanceTransform(bw, dst, cv::DIST_C, maskSizes[i], CV_32F);
        cv::normalize(dst, normDst, 0, 1, cv::NORM_MINMAX);
        dsts.push_back(dst);
        normDsts.push_back(normDst);
    }
    
    for (int i = 0; i < normDsts.size(); i++) {
        cv::imshow(winNames[i], normDsts[i]);
    }
    
    std::vector<std::string> marks = {"3 vs 5", "5 vs PRECISE", "3 vs PRECISE"};
    std::vector<cv::Mat> diffs;
    diffs.push_back(dsts[0] - dsts[1]);
    diffs.push_back(dsts[1] - dsts[2]);
    diffs.push_back(dsts[0] - dsts[2]);
    for (int i = 0; i < diffs.size(); i++) {
        std::vector<cv::Point> nonzero;
        cv::findNonZero(diffs[i], nonzero);
        std::cout << marks[i] << ": " << nonzero.size() << std::endl;
    }
    */
    /*
    cv::Mat dstL18U;
    cv::distanceTransform(bw, dstL18U, cv::DIST_L1, cv::DIST_MASK_3, CV_8U);
    bool hasDiff = false;
    for (int i = 0; i < dstL18U.rows; i++) {
        for (int j = 0; j < dstL18U.cols; j++) {
            uchar uc = dstL18U.ptr<uchar>(i)[j];
            float f = dsts[0].ptr<float>(i)[j];
            float ff = f - uc;
            if (fabs(ff) > FLT_EPSILON) {
                hasDiff = true;
                break;
            }
        }
    }
    
    std::cout << (hasDiff ? "different" : "same") << std::endl;
    */
    cv::waitKey(0);
    
    
    return 0;
}
